### Copy thy neighbour: Spatial interdependences in the democracy-repression nexus ###
### Journal of Human Rights ###
### Created by Roman-Gabriel Olar on 01/05/2019###
### Last modified on 25/08/2023###

### Load packages
require(foreign)
require(sandwich)
require(lmtest)
require(MASS)
require(Hmisc)
require(stargazer)
require(sciplot)
require(lme4)

### Clear space
rm(list = ls())


### Define function to combine multiple model estimates (this incorporates uncertainty from the latent variables)
milm <- function(fml, midata){
  xx <- terms(as.formula(fml))
  lms <- matrix(data=NA, nrow=(length(attr(xx, "term.labels")) + 1), ncol=length(midata))
  ses <- matrix(data=NA, nrow=(length(attr(xx, "term.labels")) + 1), ncol=length(midata))
  vcovs <- list()
  for(i in 1:length(midata)){
    tmp <- lm(formula=as.formula(fml), data=midata[[i]])
    lms[,i] <- tmp$coefficients
    ses[,i] <- sqrt(diag(vcov(tmp)))
    vcovs[[i]] <- vcov(tmp)
  }
  par.est <- apply(lms, 1, mean)
  se.within <- apply(ses, 1, mean)
  se.between <- apply(lms, 1, var)
  se.est <- sqrt(se.within^2 + se.between*(1 + (1/length(midata))))
  list("terms"=names(tmp$coefficients), "beta" = par.est, "SE"=se.est, "vcovs"=vcovs,"coefs" = lms)    
}

### Load data. 
setwd("set your own working directory")

dat <- read.dta("MainData.dta")

attach(dat)
datfh <- na.omit(data.frame(theta_mean, theta_sd, l_fariss, l_fariss_sd, SLreg_fariss_t1, SLreg_fariss_sd, v2cltort_sd, v2clkill_sd,
                            l_polity, l_polyarchy, l_polyarchy_sd, SLreg_polyarchy_t1, SLreg_polyarchy_sd, l_uds, l_boix, l_cgv, l_gwf_democracy, 
                            SLreg_polity, SLreg_gwf, SLreg_polyarchy_t1, SLreg_boix_t1, SLreg_cgv_t1,
                            SLreg_gdp_gro, SLreg_gdppcln, SLreg_lnpop, SLreg_civilwar, SLreg_intwar,
                            SL_def_fariss, SL_def_fariss_sd, SL_dem_fariss, SL_dem_fariss_sd,
                            SLreg_clkill_sd, SLreg_clkill_t1, SLreg_cltort_sd, SLreg_cltort_t1, l_v2clkill_sd, l_v2clkill, l_cltort_sd, l_cltort, v2cltort, v2clkill,
                            l_uds2, l_uds_sd, SLreg_uds_t1, SLreg_uds_sd, l_gdp_gro, l_lngdppc, l_lnpop, civil_war, international_war, l_resdep,
                            ccode_fe1, ccode_fe2, ccode_fe3, ccode_fe4, ccode_fe5, ccode_fe6, ccode_fe7, ccode_fe8, ccode_fe9, 
                            ccode_fe10, ccode_fe11, ccode_fe12, ccode_fe13, ccode_fe14, ccode_fe15, ccode_fe16, ccode_fe17, 
                            ccode_fe18, ccode_fe19, ccode_fe20, ccode_fe21, ccode_fe22, ccode_fe23, ccode_fe24, ccode_fe25, 
                            ccode_fe26, ccode_fe27, ccode_fe28, ccode_fe29, ccode_fe30, ccode_fe31, ccode_fe32, ccode_fe33, 
                            ccode_fe34, ccode_fe35, ccode_fe36, ccode_fe37, ccode_fe38, ccode_fe39, ccode_fe40, ccode_fe41, 
                            ccode_fe42, ccode_fe43, ccode_fe44, ccode_fe45, ccode_fe46, ccode_fe47, ccode_fe48, ccode_fe49, 
                            ccode_fe50, ccode_fe51, ccode_fe52, ccode_fe53, ccode_fe54, ccode_fe55, ccode_fe56, ccode_fe57,
                            ccode_fe58, ccode_fe59, ccode_fe60, ccode_fe61, ccode_fe62, ccode_fe63, ccode_fe64, ccode_fe65, 
                            ccode_fe66, ccode_fe67, ccode_fe68, ccode_fe69, ccode_fe70, ccode_fe71, ccode_fe72, ccode_fe73, 
                            ccode_fe74, ccode_fe75, ccode_fe76, ccode_fe77, ccode_fe78, ccode_fe79, ccode_fe80, ccode_fe81, 
                            ccode_fe82, ccode_fe83, ccode_fe84, ccode_fe85, ccode_fe86, ccode_fe87, ccode_fe88, ccode_fe89, 
                            ccode_fe90, ccode_fe91, ccode_fe92, ccode_fe93, ccode_fe94, ccode_fe95, ccode_fe96, ccode_fe97, 
                            ccode_fe98, ccode_fe99, ccode_fe100, ccode_fe101, ccode_fe102, ccode_fe103, ccode_fe104, ccode_fe105, 
                            ccode_fe106, ccode_fe107, ccode_fe108, ccode_fe109, ccode_fe110, ccode_fe111, ccode_fe112, ccode_fe113, 
                            ccode_fe114, ccode_fe115, ccode_fe116, ccode_fe117, ccode_fe118, ccode_fe119, ccode_fe120, ccode_fe121, 
                            ccode_fe122, ccode_fe123, ccode_fe124, ccode_fe125, ccode_fe126, ccode_fe127, ccode_fe128, ccode_fe129, 
                            ccode_fe130, ccode_fe131, ccode_fe132, ccode_fe133, ccode_fe134, ccode_fe135, ccode_fe136, ccode_fe137, 
                            ccode_fe138, y_fe1, y_fe2, y_fe3, y_fe4, y_fe5, y_fe6, 
                            y_fe7, y_fe8, y_fe9, y_fe10, y_fe11, y_fe12, y_fe13, y_fe14, y_fe15, y_fe16, y_fe17, y_fe18, y_fe19, y_fe20, 
                            y_fe21, y_fe22, y_fe23, y_fe24, y_fe25, y_fe26, y_fe27, y_fe28, y_fe29, y_fe30, y_fe31, y_fe32, y_fe33, y_fe34, 
                            y_fe35, y_fe36, y_fe37, y_fe38, y_fe39, y_fe40, y_fe41, y_fe42, y_fe43, y_fe44, y_fe45, y_fe46, y_fe47, y_fe48, 
                            y_fe49, y_fe50, y_fe51, y_fe52, y_fe53, y_fe54, y_fe55, y_fe56, y_fe57, y_fe58, y_fe59, y_fe60, y_fe61))

dim(datfh)


### Setup data to run models where we account for uncertainty
set.seed(32261991)
newdata <- list()
for(i in 1:1000){
  newdata[[i]] <- datfh
  newdata[[i]]$draw <- rnorm(cbind(rep(1,nrow(datfh))), mean=datfh$theta_mean, sd=datfh$theta_sd)
  newdata[[i]]$draw.lag <- rnorm(cbind(rep(1,nrow(datfh))), mean=datfh$l_fariss, sd=datfh$l_fariss_sd)
  newdata[[i]]$draw.lag.SLfariss <- rnorm(cbind(rep(1,nrow(datfh))), mean=datfh$SLreg_fariss_t1, sd=datfh$SLreg_fariss_sd)
  newdata[[i]]$draw.lag.SLuds <- rnorm(cbind(rep(1,nrow(datfh))), mean=datfh$SLreg_uds_t1, sd=datfh$SLreg_uds_sd)
  newdata[[i]]$draw.lag.uds <- rnorm(cbind(rep(1,nrow(datfh))), mean=datfh$l_uds2, sd=datfh$l_uds_sd)
  newdata[[i]]$draw.clkill <- rnorm(cbind(rep(1,nrow(datfh))), mean=datfh$v2clkill, sd=datfh$v2clkill_sd)
  newdata[[i]]$draw.lag.clkill <- rnorm(cbind(rep(1,nrow(datfh))), mean=datfh$l_v2clkill, sd=datfh$l_v2clkill_sd)
  newdata[[i]]$draw.cltort <- rnorm(cbind(rep(1,nrow(datfh))), mean=datfh$v2cltort, sd=datfh$v2cltort_sd)
  newdata[[i]]$draw.lag.cltort <- rnorm(cbind(rep(1,nrow(datfh))), mean=datfh$l_cltort, sd=datfh$l_cltort_sd)
  newdata[[i]]$draw.lag.SLkill <- rnorm(cbind(rep(1,nrow(datfh))), mean=datfh$SLreg_clkill_t1, sd=datfh$SLreg_clkill_sd)
  newdata[[i]]$draw.lag.SLtort <- rnorm(cbind(rep(1,nrow(datfh))), mean=datfh$SLreg_cltort_t1, sd=datfh$SLreg_cltort_sd)
  newdata[[i]]$draw.lag.SLpolyarchy <- rnorm(cbind(rep(1,nrow(datfh))), mean=datfh$SLreg_polyarchy_t1, sd=datfh$SLreg_polyarchy_sd)
  newdata[[i]]$draw.lag.polyarchy <- rnorm(cbind(rep(1,nrow(datfh))), mean=datfh$l_polyarchy, sd=datfh$l_polyarchy_sd)
  newdata[[i]]$draw.lag.SLdem_fariss <- rnorm(cbind(rep(1,nrow(datfh))), mean=datfh$SL_dem_fariss, sd=datfh$SL_dem_fariss_sd)
  newdata[[i]]$draw.lag.SLdef_fariss <- rnorm(cbind(rep(1,nrow(datfh))), mean=datfh$SL_def_fariss, sd=datfh$SL_def_fariss_sd)
  }

### Tables A3 and A4 ###

##Model 1 UDS - Non-spatial (NS)
FML_NS <-"draw ~ draw.lag + draw.lag.uds + 
l_gdp_gro + l_lngdppc + l_lnpop + civil_war + international_war + l_resdep"

model_ns <- milm(fml=FML_NS, midata=newdata)
df <- 5870

# create object for printing during the loop and then print the results to screen
model1<- data.frame(model_ns$terms, model_ns$beta, model_ns$SE, model_ns$beta/model_ns$SE)
colnames(model1) <- c("Variable","Est. Coef.","SE","t-statistic")

print(paste("Model Output"))
print(model1)


#Model 3 - Spatial Autoregressive (SAR)
FML_SAR <-"draw ~ draw.lag + draw.lag.SLfariss + draw.lag.uds + 
l_gdp_gro + l_lngdppc + l_lnpop + civil_war + international_war + l_resdep + ccode_fe1+ 
ccode_fe2+ ccode_fe3+ ccode_fe4 + ccode_fe5+ ccode_fe6+ ccode_fe7+ ccode_fe8+ ccode_fe9+ 
ccode_fe10+ ccode_fe11+ ccode_fe12+ ccode_fe13+ ccode_fe14+ ccode_fe15+ ccode_fe16+ ccode_fe17+ 
ccode_fe18+ ccode_fe19+ ccode_fe20+ ccode_fe21+ ccode_fe22+ ccode_fe23+ ccode_fe24+ ccode_fe25+ 
ccode_fe26+ ccode_fe27+ ccode_fe28+ ccode_fe29+ ccode_fe30+ ccode_fe31+ ccode_fe32+ ccode_fe33+ 
ccode_fe34+ ccode_fe35+ ccode_fe36+ ccode_fe37+ ccode_fe38+ ccode_fe39+ ccode_fe40+ ccode_fe41+ 
ccode_fe42+ ccode_fe43+ ccode_fe44+ ccode_fe45+ ccode_fe46+ ccode_fe47+ ccode_fe48+ ccode_fe49+ 
ccode_fe50+ ccode_fe51+ ccode_fe52+ ccode_fe53+ ccode_fe54+ ccode_fe55+ ccode_fe56+ ccode_fe57+
ccode_fe58+ ccode_fe59+ ccode_fe60+ ccode_fe61+ ccode_fe62+ ccode_fe63+ ccode_fe64+ ccode_fe65+ 
ccode_fe66+ ccode_fe67+ ccode_fe68+ ccode_fe69+ ccode_fe70+ ccode_fe71+ ccode_fe72+ ccode_fe73+ 
ccode_fe74+ ccode_fe75+ ccode_fe76+ ccode_fe77+ ccode_fe78+ ccode_fe79+ ccode_fe80+ ccode_fe81+ 
ccode_fe82+ ccode_fe83+ ccode_fe85+ ccode_fe86+ ccode_fe87+ ccode_fe88+ ccode_fe89+ 
ccode_fe90+ ccode_fe91+ ccode_fe92+ ccode_fe93+ ccode_fe94+ ccode_fe95+ ccode_fe96+ ccode_fe97+ 
ccode_fe98+ ccode_fe99+ ccode_fe100+ ccode_fe101+ ccode_fe102+ ccode_fe103+ ccode_fe104+ ccode_fe105+ 
ccode_fe106+ ccode_fe107+ ccode_fe108+ ccode_fe109+ ccode_fe110+ ccode_fe111+ ccode_fe112+ ccode_fe113+ 
ccode_fe114+ ccode_fe115+ ccode_fe116+ ccode_fe117+ ccode_fe118+ ccode_fe119+ ccode_fe120+ ccode_fe121+ 
ccode_fe122+ ccode_fe123+ ccode_fe124+ ccode_fe125+ ccode_fe126+ ccode_fe127+ ccode_fe128+ ccode_fe129+ 
ccode_fe130+ ccode_fe131+ ccode_fe132+ ccode_fe133+ ccode_fe134+ ccode_fe135+ ccode_fe136+ ccode_fe137+ 
ccode_fe138+ y_fe1+ y_fe2+ y_fe3+ y_fe5+ y_fe6+ 
y_fe7+ y_fe8+ y_fe9+ y_fe10+ y_fe11+ y_fe12+ y_fe13+ y_fe14+ y_fe15+ y_fe16+ y_fe17+ y_fe18+ y_fe19+ y_fe20+ 
y_fe21+ y_fe22+ y_fe23+ y_fe24+ y_fe25+ y_fe26+ y_fe27+ y_fe28+ y_fe29+ y_fe30+ y_fe31+ y_fe32+ y_fe33+ y_fe34+ 
y_fe35+ y_fe36+ y_fe37+ y_fe38+ y_fe39+ y_fe40+ y_fe41+ y_fe42+ y_fe43+ y_fe44+ y_fe45+ y_fe46+ y_fe47+ y_fe48+ 
y_fe49+ y_fe50+ y_fe51+ y_fe52+ y_fe53+ y_fe54+ y_fe55+ y_fe56+ y_fe57+ y_fe58+ y_fe59+ y_fe60+ y_fe61"

model_sar <- milm(fml=FML_SAR, midata=newdata)

model3<- data.frame(model_sar$terms, model_sar$beta, model_sar$SE, model_sar$beta/model_sar$SE)
colnames(model3) <- c("Variable","Est. Coef.","SE","z-score")
print(paste("Model Output"))
print(model3)


#Model 4 - Spatial Lag of X (SLX)

FML_SLX <-"draw ~ draw.lag + draw.lag.SLuds + draw.lag.uds + 
l_gdp_gro + l_lngdppc + l_lnpop + civil_war + international_war + l_resdep + ccode_fe1+ 
ccode_fe2+ ccode_fe3+ ccode_fe4 + ccode_fe5+ ccode_fe6+ ccode_fe7+ ccode_fe8+ ccode_fe9+ 
ccode_fe10+ ccode_fe11+ ccode_fe12+ ccode_fe13+ ccode_fe14+ ccode_fe15+ ccode_fe16+ ccode_fe17+ 
ccode_fe18+ ccode_fe19+ ccode_fe20+ ccode_fe21+ ccode_fe22+ ccode_fe23+ ccode_fe24+ ccode_fe25+ 
ccode_fe26+ ccode_fe27+ ccode_fe28+ ccode_fe29+ ccode_fe30+ ccode_fe31+ ccode_fe32+ ccode_fe33+ 
ccode_fe34+ ccode_fe35+ ccode_fe36+ ccode_fe37+ ccode_fe38+ ccode_fe39+ ccode_fe40+ ccode_fe41+ 
ccode_fe42+ ccode_fe43+ ccode_fe44+ ccode_fe45+ ccode_fe46+ ccode_fe47+ ccode_fe48+ ccode_fe49+ 
ccode_fe50+ ccode_fe51+ ccode_fe52+ ccode_fe53+ ccode_fe54+ ccode_fe55+ ccode_fe56+ ccode_fe57+
ccode_fe58+ ccode_fe59+ ccode_fe60+ ccode_fe61+ ccode_fe62+ ccode_fe63+ ccode_fe64+ ccode_fe65+ 
ccode_fe66+ ccode_fe67+ ccode_fe68+ ccode_fe69+ ccode_fe70+ ccode_fe71+ ccode_fe72+ ccode_fe73+ 
ccode_fe74+ ccode_fe75+ ccode_fe76+ ccode_fe77+ ccode_fe78+ ccode_fe79+ ccode_fe80+ ccode_fe81+ 
ccode_fe82+ ccode_fe83+ ccode_fe85+ ccode_fe86+ ccode_fe87+ ccode_fe88+ ccode_fe89+ 
ccode_fe90+ ccode_fe91+ ccode_fe92+ ccode_fe93+ ccode_fe94+ ccode_fe95+ ccode_fe96+ ccode_fe97+ 
ccode_fe98+ ccode_fe99+ ccode_fe100+ ccode_fe101+ ccode_fe102+ ccode_fe103+ ccode_fe104+ ccode_fe105+ 
ccode_fe106+ ccode_fe107+ ccode_fe108+ ccode_fe109+ ccode_fe110+ ccode_fe111+ ccode_fe112+ ccode_fe113+ 
ccode_fe114+ ccode_fe115+ ccode_fe116+ ccode_fe117+ ccode_fe118+ ccode_fe119+ ccode_fe120+ ccode_fe121+ 
ccode_fe122+ ccode_fe123+ ccode_fe124+ ccode_fe125+ ccode_fe126+ ccode_fe127+ ccode_fe128+ ccode_fe129+ 
ccode_fe130+ ccode_fe131+ ccode_fe132+ ccode_fe133+ ccode_fe134+ ccode_fe135+ ccode_fe136+ ccode_fe137+ 
ccode_fe138+ y_fe1+ y_fe2+ y_fe3+ y_fe5+ y_fe6+ 
y_fe7+ y_fe8+ y_fe9+ y_fe10+ y_fe11+ y_fe12+ y_fe13+ y_fe14+ y_fe15+ y_fe16+ y_fe17+ y_fe18+ y_fe19+ y_fe20+ 
y_fe21+ y_fe22+ y_fe23+ y_fe24+ y_fe25+ y_fe26+ y_fe27+ y_fe28+ y_fe29+ y_fe30+ y_fe31+ y_fe32+ y_fe33+ y_fe34+ 
y_fe35+ y_fe36+ y_fe37+ y_fe38+ y_fe39+ y_fe40+ y_fe41+ y_fe42+ y_fe43+ y_fe44+ y_fe45+ y_fe46+ y_fe47+ y_fe48+ 
y_fe49+ y_fe50+ y_fe51+ y_fe52+ y_fe53+ y_fe54+ y_fe55+ y_fe56+ y_fe57+ y_fe58+ y_fe59+ y_fe60+ y_fe61"


model_slx <- milm(fml=FML_SLX, midata=newdata)

model4<- data.frame(model_slx$terms, model_slx$beta, model_slx$SE, model_slx$beta/model_slx$SE)
colnames(model4) <- c("Variable","Est. Coef.","SE","z-score")
print(paste("Model Output"))
print(model4)


#Model 5 - Spatial Durbin Model (SDM)
FML_SDM <-"draw ~ draw.lag + draw.lag.SLfariss + draw.lag.SLuds + draw.lag.uds + 
l_gdp_gro + l_lngdppc + l_lnpop + civil_war + international_war + l_resdep + ccode_fe1+ 
ccode_fe2+ ccode_fe3+ ccode_fe4 + ccode_fe5+ ccode_fe6+ ccode_fe7+ ccode_fe8+ ccode_fe9+ 
ccode_fe10+ ccode_fe11+ ccode_fe12+ ccode_fe13+ ccode_fe14+ ccode_fe15+ ccode_fe16+ ccode_fe17+ 
ccode_fe18+ ccode_fe19+ ccode_fe20+ ccode_fe21+ ccode_fe22+ ccode_fe23+ ccode_fe24+ ccode_fe25+ 
ccode_fe26+ ccode_fe27+ ccode_fe28+ ccode_fe29+ ccode_fe30+ ccode_fe31+ ccode_fe32+ ccode_fe33+ 
ccode_fe34+ ccode_fe35+ ccode_fe36+ ccode_fe37+ ccode_fe38+ ccode_fe39+ ccode_fe40+ ccode_fe41+ 
ccode_fe42+ ccode_fe43+ ccode_fe44+ ccode_fe45+ ccode_fe46+ ccode_fe47+ ccode_fe48+ ccode_fe49+ 
ccode_fe50+ ccode_fe51+ ccode_fe52+ ccode_fe53+ ccode_fe54+ ccode_fe55+ ccode_fe56+ ccode_fe57+
ccode_fe58+ ccode_fe59+ ccode_fe60+ ccode_fe61+ ccode_fe62+ ccode_fe63+ ccode_fe64+ ccode_fe65+ 
ccode_fe66+ ccode_fe67+ ccode_fe68+ ccode_fe69+ ccode_fe70+ ccode_fe71+ ccode_fe72+ ccode_fe73+ 
ccode_fe74+ ccode_fe75+ ccode_fe76+ ccode_fe77+ ccode_fe78+ ccode_fe79+ ccode_fe80+ ccode_fe81+ 
ccode_fe82+ ccode_fe83+ ccode_fe85+ ccode_fe86+ ccode_fe87+ ccode_fe88+ ccode_fe89+ 
ccode_fe90+ ccode_fe91+ ccode_fe92+ ccode_fe93+ ccode_fe94+ ccode_fe95+ ccode_fe96+ ccode_fe97+ 
ccode_fe98+ ccode_fe99+ ccode_fe100+ ccode_fe101+ ccode_fe102+ ccode_fe103+ ccode_fe104+ ccode_fe105+ 
ccode_fe106+ ccode_fe107+ ccode_fe108+ ccode_fe109+ ccode_fe110+ ccode_fe111+ ccode_fe112+ ccode_fe113+ 
ccode_fe114+ ccode_fe115+ ccode_fe116+ ccode_fe117+ ccode_fe118+ ccode_fe119+ ccode_fe120+ ccode_fe121+ 
ccode_fe122+ ccode_fe123+ ccode_fe124+ ccode_fe125+ ccode_fe126+ ccode_fe127+ ccode_fe128+ ccode_fe129+ 
ccode_fe130+ ccode_fe131+ ccode_fe132+ ccode_fe133+ ccode_fe134+ ccode_fe135+ ccode_fe136+ ccode_fe137+ 
ccode_fe138+ y_fe1+ y_fe2+ y_fe3+ y_fe5+ y_fe6+ 
y_fe7+ y_fe8+ y_fe9+ y_fe10+ y_fe11+ y_fe12+ y_fe13+ y_fe14+ y_fe15+ y_fe16+ y_fe17+ y_fe18+ y_fe19+ y_fe20+ 
y_fe21+ y_fe22+ y_fe23+ y_fe24+ y_fe25+ y_fe26+ y_fe27+ y_fe28+ y_fe29+ y_fe30+ y_fe31+ y_fe32+ y_fe33+ y_fe34+ 
y_fe35+ y_fe36+ y_fe37+ y_fe38+ y_fe39+ y_fe40+ y_fe41+ y_fe42+ y_fe43+ y_fe44+ y_fe45+ y_fe46+ y_fe47+ y_fe48+ 
y_fe49+ y_fe50+ y_fe51+ y_fe52+ y_fe53+ y_fe54+ y_fe55+ y_fe56+ y_fe57+ y_fe58+ y_fe59+ y_fe60+ y_fe61"


model_sdm <- milm(fml=FML_SDM, midata=newdata)

model5<- data.frame(model_sdm$terms, model_sdm$beta, model_sdm$SE, model_sdm$beta/model_sdm$SE)
colnames(model5) <- c("Variable","Est. Coef.","SE","z-score")
print(paste("Model Output"))
print(model5)


### Tables A5 and A6 ###


##Model Polity - Non-spatial (NS)
FML_NS_POL <-"draw ~ draw.lag + l_polity +
l_gdp_gro + l_lngdppc + l_lnpop + civil_war + international_war + l_resdep"

model_ns_pol <- milm(fml=FML_NS_POL, midata=newdata)
df <- 5870

# create object for printing during the loop and then print the results to screen
model_pol<- data.frame(model_ns_pol$terms, model_ns_pol$beta, model_ns_pol$SE, model_ns_pol$beta/model_ns_pol$SE)
colnames(model_pol) <- c("Variable","Est. Coef.","SE","z-score")

print(paste("Polity NS"))
print(model_pol)

#Model Polity - SAR)
FML_SAR_POL <-"draw ~ draw.lag + draw.lag.SLfariss + l_polity +
l_gdp_gro + l_lngdppc + l_lnpop + civil_war + international_war + l_resdep + ccode_fe1+ 
ccode_fe2+ ccode_fe3+ ccode_fe4 + ccode_fe5+ ccode_fe6+ ccode_fe7+ ccode_fe8+ ccode_fe9+ 
ccode_fe10+ ccode_fe11+ ccode_fe12+ ccode_fe13+ ccode_fe14+ ccode_fe15+ ccode_fe16+ ccode_fe17+ 
ccode_fe18+ ccode_fe19+ ccode_fe20+ ccode_fe21+ ccode_fe22+ ccode_fe23+ ccode_fe24+ ccode_fe25+ 
ccode_fe26+ ccode_fe27+ ccode_fe28+ ccode_fe29+ ccode_fe30+ ccode_fe31+ ccode_fe32+ ccode_fe33+ 
ccode_fe34+ ccode_fe35+ ccode_fe36+ ccode_fe37+ ccode_fe38+ ccode_fe39+ ccode_fe40+ ccode_fe41+ 
ccode_fe42+ ccode_fe43+ ccode_fe44+ ccode_fe45+ ccode_fe46+ ccode_fe47+ ccode_fe48+ ccode_fe49+ 
ccode_fe50+ ccode_fe51+ ccode_fe52+ ccode_fe53+ ccode_fe54+ ccode_fe55+ ccode_fe56+ ccode_fe57+
ccode_fe58+ ccode_fe59+ ccode_fe60+ ccode_fe61+ ccode_fe62+ ccode_fe63+ ccode_fe64+ ccode_fe65+ 
ccode_fe66+ ccode_fe67+ ccode_fe68+ ccode_fe69+ ccode_fe70+ ccode_fe71+ ccode_fe72+ ccode_fe73+ 
ccode_fe74+ ccode_fe75+ ccode_fe76+ ccode_fe77+ ccode_fe78+ ccode_fe79+ ccode_fe80+ ccode_fe81+ 
ccode_fe82+ ccode_fe83+ ccode_fe85+ ccode_fe86+ ccode_fe87+ ccode_fe88+ ccode_fe89+ 
ccode_fe90+ ccode_fe91+ ccode_fe92+ ccode_fe93+ ccode_fe94+ ccode_fe95+ ccode_fe96+ ccode_fe97+ 
ccode_fe98+ ccode_fe99+ ccode_fe100+ ccode_fe101+ ccode_fe102+ ccode_fe103+ ccode_fe104+ ccode_fe105+ 
ccode_fe106+ ccode_fe107+ ccode_fe108+ ccode_fe109+ ccode_fe110+ ccode_fe111+ ccode_fe112+ ccode_fe113+ 
ccode_fe114+ ccode_fe115+ ccode_fe116+ ccode_fe117+ ccode_fe118+ ccode_fe119+ ccode_fe120+ ccode_fe121+ 
ccode_fe122+ ccode_fe123+ ccode_fe124+ ccode_fe125+ ccode_fe126+ ccode_fe127+ ccode_fe128+ ccode_fe129+ 
ccode_fe130+ ccode_fe131+ ccode_fe132+ ccode_fe133+ ccode_fe134+ ccode_fe135+ ccode_fe136+ ccode_fe137+ 
ccode_fe138+ y_fe1+ y_fe2+ y_fe3+ y_fe5+ y_fe6+ 
y_fe7+ y_fe8+ y_fe9+ y_fe10+ y_fe11+ y_fe12+ y_fe13+ y_fe14+ y_fe15+ y_fe16+ y_fe17+ y_fe18+ y_fe19+ y_fe20+ 
y_fe21+ y_fe22+ y_fe23+ y_fe24+ y_fe25+ y_fe26+ y_fe27+ y_fe28+ y_fe29+ y_fe30+ y_fe31+ y_fe32+ y_fe33+ y_fe34+ 
y_fe35+ y_fe36+ y_fe37+ y_fe38+ y_fe39+ y_fe40+ y_fe41+ y_fe42+ y_fe43+ y_fe44+ y_fe45+ y_fe46+ y_fe47+ y_fe48+ 
y_fe49+ y_fe50+ y_fe51+ y_fe52+ y_fe53+ y_fe54+ y_fe55+ y_fe56+ y_fe57+ y_fe58+ y_fe59+ y_fe60+ y_fe61"

model_sar_pol <- milm(fml=FML_SAR_POL, midata=newdata)
model_sar_pol<- data.frame(model_sar_pol$terms, model_sar_pol$beta, model_sar_pol$SE, model_sar_pol$beta/model_sar_pol$SE)
colnames(model_sar_pol) <- c("Variable","Est. Coef.","SE","z-score")

print(paste("POlity SAR"))
print(model_sar_pol)


#Model Polity - SLX
FML_SLX_POL <-"draw ~ draw.lag + SLreg_polity + l_polity +
l_gdp_gro + l_lngdppc + l_lnpop + civil_war + international_war + l_resdep + ccode_fe1+ 
ccode_fe2+ ccode_fe3+ ccode_fe4 + ccode_fe5+ ccode_fe6+ ccode_fe7+ ccode_fe8+ ccode_fe9+ 
ccode_fe10+ ccode_fe11+ ccode_fe12+ ccode_fe13+ ccode_fe14+ ccode_fe15+ ccode_fe16+ ccode_fe17+ 
ccode_fe18+ ccode_fe19+ ccode_fe20+ ccode_fe21+ ccode_fe22+ ccode_fe23+ ccode_fe24+ ccode_fe25+ 
ccode_fe26+ ccode_fe27+ ccode_fe28+ ccode_fe29+ ccode_fe30+ ccode_fe31+ ccode_fe32+ ccode_fe33+ 
ccode_fe34+ ccode_fe35+ ccode_fe36+ ccode_fe37+ ccode_fe38+ ccode_fe39+ ccode_fe40+ ccode_fe41+ 
ccode_fe42+ ccode_fe43+ ccode_fe44+ ccode_fe45+ ccode_fe46+ ccode_fe47+ ccode_fe48+ ccode_fe49+ 
ccode_fe50+ ccode_fe51+ ccode_fe52+ ccode_fe53+ ccode_fe54+ ccode_fe55+ ccode_fe56+ ccode_fe57+
ccode_fe58+ ccode_fe59+ ccode_fe60+ ccode_fe61+ ccode_fe62+ ccode_fe63+ ccode_fe64+ ccode_fe65+ 
ccode_fe66+ ccode_fe67+ ccode_fe68+ ccode_fe69+ ccode_fe70+ ccode_fe71+ ccode_fe72+ ccode_fe73+ 
ccode_fe74+ ccode_fe75+ ccode_fe76+ ccode_fe77+ ccode_fe78+ ccode_fe79+ ccode_fe80+ ccode_fe81+ 
ccode_fe82+ ccode_fe83+ ccode_fe85+ ccode_fe86+ ccode_fe87+ ccode_fe88+ ccode_fe89+ 
ccode_fe90+ ccode_fe91+ ccode_fe92+ ccode_fe93+ ccode_fe94+ ccode_fe95+ ccode_fe96+ ccode_fe97+ 
ccode_fe98+ ccode_fe99+ ccode_fe100+ ccode_fe101+ ccode_fe102+ ccode_fe103+ ccode_fe104+ ccode_fe105+ 
ccode_fe106+ ccode_fe107+ ccode_fe108+ ccode_fe109+ ccode_fe110+ ccode_fe111+ ccode_fe112+ ccode_fe113+ 
ccode_fe114+ ccode_fe115+ ccode_fe116+ ccode_fe117+ ccode_fe118+ ccode_fe119+ ccode_fe120+ ccode_fe121+ 
ccode_fe122+ ccode_fe123+ ccode_fe124+ ccode_fe125+ ccode_fe126+ ccode_fe127+ ccode_fe128+ ccode_fe129+ 
ccode_fe130+ ccode_fe131+ ccode_fe132+ ccode_fe133+ ccode_fe134+ ccode_fe135+ ccode_fe136+ ccode_fe137+ 
ccode_fe138+ y_fe1+ y_fe2+ y_fe3+ y_fe5+ y_fe6+ 
y_fe7+ y_fe8+ y_fe9+ y_fe10+ y_fe11+ y_fe12+ y_fe13+ y_fe14+ y_fe15+ y_fe16+ y_fe17+ y_fe18+ y_fe19+ y_fe20+ 
y_fe21+ y_fe22+ y_fe23+ y_fe24+ y_fe25+ y_fe26+ y_fe27+ y_fe28+ y_fe29+ y_fe30+ y_fe31+ y_fe32+ y_fe33+ y_fe34+ 
y_fe35+ y_fe36+ y_fe37+ y_fe38+ y_fe39+ y_fe40+ y_fe41+ y_fe42+ y_fe43+ y_fe44+ y_fe45+ y_fe46+ y_fe47+ y_fe48+ 
y_fe49+ y_fe50+ y_fe51+ y_fe52+ y_fe53+ y_fe54+ y_fe55+ y_fe56+ y_fe57+ y_fe58+ y_fe59+ y_fe60+ y_fe61"

model_slx_pol <- milm(fml=FML_SLX_POL, midata=newdata)
model_slx_pol<- data.frame(model_slx_pol$terms, model_slx_pol$beta, model_slx_pol$SE, model_slx_pol$beta/model_slx_pol$SE)
colnames(model_slx_pol) <- c("Variable","Est. Coef.","SE","z-score")

print(paste("POlity SLX"))
print(model_slx_pol)



#Model Polity - Spatial Durbin Model (SDM)
FML_SDM_POL <-"draw ~ draw.lag + draw.lag.SLfariss + SLreg_polity + l_polity +
l_gdp_gro + l_lngdppc + l_lnpop + civil_war + international_war + l_resdep + ccode_fe1+ 
ccode_fe2+ ccode_fe3+ ccode_fe4 + ccode_fe5+ ccode_fe6+ ccode_fe7+ ccode_fe8+ ccode_fe9+ 
ccode_fe10+ ccode_fe11+ ccode_fe12+ ccode_fe13+ ccode_fe14+ ccode_fe15+ ccode_fe16+ ccode_fe17+ 
ccode_fe18+ ccode_fe19+ ccode_fe20+ ccode_fe21+ ccode_fe22+ ccode_fe23+ ccode_fe24+ ccode_fe25+ 
ccode_fe26+ ccode_fe27+ ccode_fe28+ ccode_fe29+ ccode_fe30+ ccode_fe31+ ccode_fe32+ ccode_fe33+ 
ccode_fe34+ ccode_fe35+ ccode_fe36+ ccode_fe37+ ccode_fe38+ ccode_fe39+ ccode_fe40+ ccode_fe41+ 
ccode_fe42+ ccode_fe43+ ccode_fe44+ ccode_fe45+ ccode_fe46+ ccode_fe47+ ccode_fe48+ ccode_fe49+ 
ccode_fe50+ ccode_fe51+ ccode_fe52+ ccode_fe53+ ccode_fe54+ ccode_fe55+ ccode_fe56+ ccode_fe57+
ccode_fe58+ ccode_fe59+ ccode_fe60+ ccode_fe61+ ccode_fe62+ ccode_fe63+ ccode_fe64+ ccode_fe65+ 
ccode_fe66+ ccode_fe67+ ccode_fe68+ ccode_fe69+ ccode_fe70+ ccode_fe71+ ccode_fe72+ ccode_fe73+ 
ccode_fe74+ ccode_fe75+ ccode_fe76+ ccode_fe77+ ccode_fe78+ ccode_fe79+ ccode_fe80+ ccode_fe81+ 
ccode_fe82+ ccode_fe83+ ccode_fe85+ ccode_fe86+ ccode_fe87+ ccode_fe88+ ccode_fe89+ 
ccode_fe90+ ccode_fe91+ ccode_fe92+ ccode_fe93+ ccode_fe94+ ccode_fe95+ ccode_fe96+ ccode_fe97+ 
ccode_fe98+ ccode_fe99+ ccode_fe100+ ccode_fe101+ ccode_fe102+ ccode_fe103+ ccode_fe104+ ccode_fe105+ 
ccode_fe106+ ccode_fe107+ ccode_fe108+ ccode_fe109+ ccode_fe110+ ccode_fe111+ ccode_fe112+ ccode_fe113+ 
ccode_fe114+ ccode_fe115+ ccode_fe116+ ccode_fe117+ ccode_fe118+ ccode_fe119+ ccode_fe120+ ccode_fe121+ 
ccode_fe122+ ccode_fe123+ ccode_fe124+ ccode_fe125+ ccode_fe126+ ccode_fe127+ ccode_fe128+ ccode_fe129+ 
ccode_fe130+ ccode_fe131+ ccode_fe132+ ccode_fe133+ ccode_fe134+ ccode_fe135+ ccode_fe136+ ccode_fe137+ 
ccode_fe138+ y_fe1+ y_fe2+ y_fe3+ y_fe5+ y_fe6+ 
y_fe7+ y_fe8+ y_fe9+ y_fe10+ y_fe11+ y_fe12+ y_fe13+ y_fe14+ y_fe15+ y_fe16+ y_fe17+ y_fe18+ y_fe19+ y_fe20+ 
y_fe21+ y_fe22+ y_fe23+ y_fe24+ y_fe25+ y_fe26+ y_fe27+ y_fe28+ y_fe29+ y_fe30+ y_fe31+ y_fe32+ y_fe33+ y_fe34+ 
y_fe35+ y_fe36+ y_fe37+ y_fe38+ y_fe39+ y_fe40+ y_fe41+ y_fe42+ y_fe43+ y_fe44+ y_fe45+ y_fe46+ y_fe47+ y_fe48+ 
y_fe49+ y_fe50+ y_fe51+ y_fe52+ y_fe53+ y_fe54+ y_fe55+ y_fe56+ y_fe57+ y_fe58+ y_fe59+ y_fe60+ y_fe61"

model_sdm_pol <- milm(fml=FML_SDM_POL, midata=newdata)
model_sdm_pol<- data.frame(model_sdm_pol$terms, model_sdm_pol$beta, model_sdm_pol$SE, model_sdm_pol$beta/model_sdm_pol$SE)
colnames(model_sdm_pol) <- c("Variable","Est. Coef.","SE","z-score")

print(paste("POlity SDM"))
print(model_sdm_pol)

### Tables A7 and A8 ###


##Model CGV - Non-spatial (NS)
FML_NS_CGV <-"draw ~ draw.lag + l_cgv +
l_gdp_gro + l_lngdppc + l_lnpop + civil_war + international_war + l_resdep"

model_ns_cgv <- milm(fml=FML_NS_CGV, midata=newdata)
df <- 5870

# create object for printing during the loop and then print the results to screen
model_cgv<- data.frame(model_ns_cgv$terms, model_ns_cgv$beta, model_ns_cgv$SE, model_ns_cgv$beta/model_ns_cgv$SE)
colnames(model_cgv) <- c("Variable","Est. Coef.","SE","z-score")

print(paste("CGV NS"))
print(model_cgv)


#Model CGV - SAR
FML_SAR_CGV <-"draw ~ draw.lag + draw.lag.SLfariss + l_cgv +
l_gdp_gro + l_lngdppc + l_lnpop + civil_war + international_war + l_resdep + ccode_fe1+ 
ccode_fe2+ ccode_fe3+ ccode_fe4 + ccode_fe5+ ccode_fe6+ ccode_fe7+ ccode_fe8+ ccode_fe9+ 
ccode_fe10+ ccode_fe11+ ccode_fe12+ ccode_fe13+ ccode_fe14+ ccode_fe15+ ccode_fe16+ ccode_fe17+ 
ccode_fe18+ ccode_fe19+ ccode_fe20+ ccode_fe21+ ccode_fe22+ ccode_fe23+ ccode_fe24+ ccode_fe25+ 
ccode_fe26+ ccode_fe27+ ccode_fe28+ ccode_fe29+ ccode_fe30+ ccode_fe31+ ccode_fe32+ ccode_fe33+ 
ccode_fe34+ ccode_fe35+ ccode_fe36+ ccode_fe37+ ccode_fe38+ ccode_fe39+ ccode_fe40+ ccode_fe41+ 
ccode_fe42+ ccode_fe43+ ccode_fe44+ ccode_fe45+ ccode_fe46+ ccode_fe47+ ccode_fe48+ ccode_fe49+ 
ccode_fe50+ ccode_fe51+ ccode_fe52+ ccode_fe53+ ccode_fe54+ ccode_fe55+ ccode_fe56+ ccode_fe57+
ccode_fe58+ ccode_fe59+ ccode_fe60+ ccode_fe61+ ccode_fe62+ ccode_fe63+ ccode_fe64+ ccode_fe65+ 
ccode_fe66+ ccode_fe67+ ccode_fe68+ ccode_fe69+ ccode_fe70+ ccode_fe71+ ccode_fe72+ ccode_fe73+ 
ccode_fe74+ ccode_fe75+ ccode_fe76+ ccode_fe77+ ccode_fe78+ ccode_fe79+ ccode_fe80+ ccode_fe81+ 
ccode_fe82+ ccode_fe83+ ccode_fe85+ ccode_fe86+ ccode_fe87+ ccode_fe88+ ccode_fe89+ 
ccode_fe90+ ccode_fe91+ ccode_fe92+ ccode_fe93+ ccode_fe94+ ccode_fe95+ ccode_fe96+ ccode_fe97+ 
ccode_fe98+ ccode_fe99+ ccode_fe100+ ccode_fe101+ ccode_fe102+ ccode_fe103+ ccode_fe104+ ccode_fe105+ 
ccode_fe106+ ccode_fe107+ ccode_fe108+ ccode_fe109+ ccode_fe110+ ccode_fe111+ ccode_fe112+ ccode_fe113+ 
ccode_fe114+ ccode_fe115+ ccode_fe116+ ccode_fe117+ ccode_fe118+ ccode_fe119+ ccode_fe120+ ccode_fe121+ 
ccode_fe122+ ccode_fe123+ ccode_fe124+ ccode_fe125+ ccode_fe126+ ccode_fe127+ ccode_fe128+ ccode_fe129+ 
ccode_fe130+ ccode_fe131+ ccode_fe132+ ccode_fe133+ ccode_fe134+ ccode_fe135+ ccode_fe136+ ccode_fe137+ 
ccode_fe138+ y_fe1+ y_fe2+ y_fe3+ y_fe5+ y_fe6+ 
y_fe7+ y_fe8+ y_fe9+ y_fe10+ y_fe11+ y_fe12+ y_fe13+ y_fe14+ y_fe15+ y_fe16+ y_fe17+ y_fe18+ y_fe19+ y_fe20+ 
y_fe21+ y_fe22+ y_fe23+ y_fe24+ y_fe25+ y_fe26+ y_fe27+ y_fe28+ y_fe29+ y_fe30+ y_fe31+ y_fe32+ y_fe33+ y_fe34+ 
y_fe35+ y_fe36+ y_fe37+ y_fe38+ y_fe39+ y_fe40+ y_fe41+ y_fe42+ y_fe43+ y_fe44+ y_fe45+ y_fe46+ y_fe47+ y_fe48+ 
y_fe49+ y_fe50+ y_fe51+ y_fe52+ y_fe53+ y_fe54+ y_fe55+ y_fe56+ y_fe57+ y_fe58+ y_fe59+ y_fe60+ y_fe61"

model_sar_cgv <- milm(fml=FML_SAR_CGV, midata=newdata)
model_sar_cgv<- data.frame(model_sar_cgv$terms, model_sar_cgv$beta, model_sar_cgv$SE, model_sar_cgv$beta/model_sar_cgv$SE)
colnames(model_sar_cgv) <- c("Variable","Est. Coef.","SE","z-score")

print(paste("CGV SAR"))
print(model_sar_cgv)


#Model CGV - SLX
FML_SLX_CGV <-"draw ~ draw.lag + SLreg_cgv_t1 + l_cgv +
l_gdp_gro + l_lngdppc + l_lnpop + civil_war + international_war + l_resdep + ccode_fe1+ 
ccode_fe2+ ccode_fe3+ ccode_fe4 + ccode_fe5+ ccode_fe6+ ccode_fe7+ ccode_fe8+ ccode_fe9+ 
ccode_fe10+ ccode_fe11+ ccode_fe12+ ccode_fe13+ ccode_fe14+ ccode_fe15+ ccode_fe16+ ccode_fe17+ 
ccode_fe18+ ccode_fe19+ ccode_fe20+ ccode_fe21+ ccode_fe22+ ccode_fe23+ ccode_fe24+ ccode_fe25+ 
ccode_fe26+ ccode_fe27+ ccode_fe28+ ccode_fe29+ ccode_fe30+ ccode_fe31+ ccode_fe32+ ccode_fe33+ 
ccode_fe34+ ccode_fe35+ ccode_fe36+ ccode_fe37+ ccode_fe38+ ccode_fe39+ ccode_fe40+ ccode_fe41+ 
ccode_fe42+ ccode_fe43+ ccode_fe44+ ccode_fe45+ ccode_fe46+ ccode_fe47+ ccode_fe48+ ccode_fe49+ 
ccode_fe50+ ccode_fe51+ ccode_fe52+ ccode_fe53+ ccode_fe54+ ccode_fe55+ ccode_fe56+ ccode_fe57+
ccode_fe58+ ccode_fe59+ ccode_fe60+ ccode_fe61+ ccode_fe62+ ccode_fe63+ ccode_fe64+ ccode_fe65+ 
ccode_fe66+ ccode_fe67+ ccode_fe68+ ccode_fe69+ ccode_fe70+ ccode_fe71+ ccode_fe72+ ccode_fe73+ 
ccode_fe74+ ccode_fe75+ ccode_fe76+ ccode_fe77+ ccode_fe78+ ccode_fe79+ ccode_fe80+ ccode_fe81+ 
ccode_fe82+ ccode_fe83+ ccode_fe85+ ccode_fe86+ ccode_fe87+ ccode_fe88+ ccode_fe89+ 
ccode_fe90+ ccode_fe91+ ccode_fe92+ ccode_fe93+ ccode_fe94+ ccode_fe95+ ccode_fe96+ ccode_fe97+ 
ccode_fe98+ ccode_fe99+ ccode_fe100+ ccode_fe101+ ccode_fe102+ ccode_fe103+ ccode_fe104+ ccode_fe105+ 
ccode_fe106+ ccode_fe107+ ccode_fe108+ ccode_fe109+ ccode_fe110+ ccode_fe111+ ccode_fe112+ ccode_fe113+ 
ccode_fe114+ ccode_fe115+ ccode_fe116+ ccode_fe117+ ccode_fe118+ ccode_fe119+ ccode_fe120+ ccode_fe121+ 
ccode_fe122+ ccode_fe123+ ccode_fe124+ ccode_fe125+ ccode_fe126+ ccode_fe127+ ccode_fe128+ ccode_fe129+ 
ccode_fe130+ ccode_fe131+ ccode_fe132+ ccode_fe133+ ccode_fe134+ ccode_fe135+ ccode_fe136+ ccode_fe137+ 
ccode_fe138+ y_fe1+ y_fe2+ y_fe3+ y_fe5+ y_fe6+ 
y_fe7+ y_fe8+ y_fe9+ y_fe10+ y_fe11+ y_fe12+ y_fe13+ y_fe14+ y_fe15+ y_fe16+ y_fe17+ y_fe18+ y_fe19+ y_fe20+ 
y_fe21+ y_fe22+ y_fe23+ y_fe24+ y_fe25+ y_fe26+ y_fe27+ y_fe28+ y_fe29+ y_fe30+ y_fe31+ y_fe32+ y_fe33+ y_fe34+ 
y_fe35+ y_fe36+ y_fe37+ y_fe38+ y_fe39+ y_fe40+ y_fe41+ y_fe42+ y_fe43+ y_fe44+ y_fe45+ y_fe46+ y_fe47+ y_fe48+ 
y_fe49+ y_fe50+ y_fe51+ y_fe52+ y_fe53+ y_fe54+ y_fe55+ y_fe56+ y_fe57+ y_fe58+ y_fe59+ y_fe60+ y_fe61"

model_slx_cgv <- milm(fml=FML_SLX_CGV, midata=newdata)
model_slx_cgv<- data.frame(model_slx_cgv$terms, model_slx_cgv$beta, model_slx_cgv$SE, model_slx_cgv$beta/model_slx_cgv$SE)
colnames(model_slx_cgv) <- c("Variable","Est. Coef.","SE","z-score")

print(paste("CGV SLX"))
print(model_slx_cgv)



#Model CGV - Spatial Durbin Model (SDM)
FML_SDM_CGV <-"draw ~ draw.lag + draw.lag.SLfariss + SLreg_cgv_t1 + l_cgv +
l_gdp_gro + l_lngdppc + l_lnpop + civil_war + international_war + l_resdep + ccode_fe1+ 
ccode_fe2+ ccode_fe3+ ccode_fe4 + ccode_fe5+ ccode_fe6+ ccode_fe7+ ccode_fe8+ ccode_fe9+ 
ccode_fe10+ ccode_fe11+ ccode_fe12+ ccode_fe13+ ccode_fe14+ ccode_fe15+ ccode_fe16+ ccode_fe17+ 
ccode_fe18+ ccode_fe19+ ccode_fe20+ ccode_fe21+ ccode_fe22+ ccode_fe23+ ccode_fe24+ ccode_fe25+ 
ccode_fe26+ ccode_fe27+ ccode_fe28+ ccode_fe29+ ccode_fe30+ ccode_fe31+ ccode_fe32+ ccode_fe33+ 
ccode_fe34+ ccode_fe35+ ccode_fe36+ ccode_fe37+ ccode_fe38+ ccode_fe39+ ccode_fe40+ ccode_fe41+ 
ccode_fe42+ ccode_fe43+ ccode_fe44+ ccode_fe45+ ccode_fe46+ ccode_fe47+ ccode_fe48+ ccode_fe49+ 
ccode_fe50+ ccode_fe51+ ccode_fe52+ ccode_fe53+ ccode_fe54+ ccode_fe55+ ccode_fe56+ ccode_fe57+
ccode_fe58+ ccode_fe59+ ccode_fe60+ ccode_fe61+ ccode_fe62+ ccode_fe63+ ccode_fe64+ ccode_fe65+ 
ccode_fe66+ ccode_fe67+ ccode_fe68+ ccode_fe69+ ccode_fe70+ ccode_fe71+ ccode_fe72+ ccode_fe73+ 
ccode_fe74+ ccode_fe75+ ccode_fe76+ ccode_fe77+ ccode_fe78+ ccode_fe79+ ccode_fe80+ ccode_fe81+ 
ccode_fe82+ ccode_fe83+ ccode_fe85+ ccode_fe86+ ccode_fe87+ ccode_fe88+ ccode_fe89+ 
ccode_fe90+ ccode_fe91+ ccode_fe92+ ccode_fe93+ ccode_fe94+ ccode_fe95+ ccode_fe96+ ccode_fe97+ 
ccode_fe98+ ccode_fe99+ ccode_fe100+ ccode_fe101+ ccode_fe102+ ccode_fe103+ ccode_fe104+ ccode_fe105+ 
ccode_fe106+ ccode_fe107+ ccode_fe108+ ccode_fe109+ ccode_fe110+ ccode_fe111+ ccode_fe112+ ccode_fe113+ 
ccode_fe114+ ccode_fe115+ ccode_fe116+ ccode_fe117+ ccode_fe118+ ccode_fe119+ ccode_fe120+ ccode_fe121+ 
ccode_fe122+ ccode_fe123+ ccode_fe124+ ccode_fe125+ ccode_fe126+ ccode_fe127+ ccode_fe128+ ccode_fe129+ 
ccode_fe130+ ccode_fe131+ ccode_fe132+ ccode_fe133+ ccode_fe134+ ccode_fe135+ ccode_fe136+ ccode_fe137+ 
ccode_fe138+ y_fe1+ y_fe2+ y_fe3+ y_fe5+ y_fe6+ 
y_fe7+ y_fe8+ y_fe9+ y_fe10+ y_fe11+ y_fe12+ y_fe13+ y_fe14+ y_fe15+ y_fe16+ y_fe17+ y_fe18+ y_fe19+ y_fe20+ 
y_fe21+ y_fe22+ y_fe23+ y_fe24+ y_fe25+ y_fe26+ y_fe27+ y_fe28+ y_fe29+ y_fe30+ y_fe31+ y_fe32+ y_fe33+ y_fe34+ 
y_fe35+ y_fe36+ y_fe37+ y_fe38+ y_fe39+ y_fe40+ y_fe41+ y_fe42+ y_fe43+ y_fe44+ y_fe45+ y_fe46+ y_fe47+ y_fe48+ 
y_fe49+ y_fe50+ y_fe51+ y_fe52+ y_fe53+ y_fe54+ y_fe55+ y_fe56+ y_fe57+ y_fe58+ y_fe59+ y_fe60+ y_fe61"

model_sdm_cgv <- milm(fml=FML_SDM_CGV, midata=newdata)
model_sdm_cgv<- data.frame(model_sdm_cgv$terms, model_sdm_cgv$beta, model_sdm_cgv$SE, model_sdm_cgv$beta/model_sdm_cgv$SE)
colnames(model_sdm_cgv) <- c("Variable","Est. Coef.","SE","z-score")

print(paste("CGV SDM"))
print(model_sdm_cgv)


### Tables A9 and A10 ###


##Model BOIX - Non-spatial (NS)
FML_NS_BOIX <-"draw ~ draw.lag + l_boix +
l_gdp_gro + l_lngdppc + l_lnpop + civil_war + international_war + l_resdep"

model_ns_boix <- milm(fml=FML_NS_BOIX, midata=newdata)
df <- 5870

# create object for printing during the loop and then print the results to screen
model_boix<- data.frame(model_ns_boix$terms, model_ns_boix$beta, model_ns_boix$SE, model_ns_boix$beta/model_ns_boix$SE)
colnames(model_boix) <- c("Variable","Est. Coef.","SE","z-score")

print(paste("BMR NS"))
print(model_boix)


#Model BOIX - SAR
FML_SAR_BOIX <-"draw ~ draw.lag + draw.lag.SLfariss + l_boix +
l_gdp_gro + l_lngdppc + l_lnpop + civil_war + international_war + l_resdep + ccode_fe1+ 
ccode_fe2+ ccode_fe3+ ccode_fe4 + ccode_fe5+ ccode_fe6+ ccode_fe7+ ccode_fe8+ ccode_fe9+ 
ccode_fe10+ ccode_fe11+ ccode_fe12+ ccode_fe13+ ccode_fe14+ ccode_fe15+ ccode_fe16+ ccode_fe17+ 
ccode_fe18+ ccode_fe19+ ccode_fe20+ ccode_fe21+ ccode_fe22+ ccode_fe23+ ccode_fe24+ ccode_fe25+ 
ccode_fe26+ ccode_fe27+ ccode_fe28+ ccode_fe29+ ccode_fe30+ ccode_fe31+ ccode_fe32+ ccode_fe33+ 
ccode_fe34+ ccode_fe35+ ccode_fe36+ ccode_fe37+ ccode_fe38+ ccode_fe39+ ccode_fe40+ ccode_fe41+ 
ccode_fe42+ ccode_fe43+ ccode_fe44+ ccode_fe45+ ccode_fe46+ ccode_fe47+ ccode_fe48+ ccode_fe49+ 
ccode_fe50+ ccode_fe51+ ccode_fe52+ ccode_fe53+ ccode_fe54+ ccode_fe55+ ccode_fe56+ ccode_fe57+
ccode_fe58+ ccode_fe59+ ccode_fe60+ ccode_fe61+ ccode_fe62+ ccode_fe63+ ccode_fe64+ ccode_fe65+ 
ccode_fe66+ ccode_fe67+ ccode_fe68+ ccode_fe69+ ccode_fe70+ ccode_fe71+ ccode_fe72+ ccode_fe73+ 
ccode_fe74+ ccode_fe75+ ccode_fe76+ ccode_fe77+ ccode_fe78+ ccode_fe79+ ccode_fe80+ ccode_fe81+ 
ccode_fe82+ ccode_fe83+ ccode_fe85+ ccode_fe86+ ccode_fe87+ ccode_fe88+ ccode_fe89+ 
ccode_fe90+ ccode_fe91+ ccode_fe92+ ccode_fe93+ ccode_fe94+ ccode_fe95+ ccode_fe96+ ccode_fe97+ 
ccode_fe98+ ccode_fe99+ ccode_fe100+ ccode_fe101+ ccode_fe102+ ccode_fe103+ ccode_fe104+ ccode_fe105+ 
ccode_fe106+ ccode_fe107+ ccode_fe108+ ccode_fe109+ ccode_fe110+ ccode_fe111+ ccode_fe112+ ccode_fe113+ 
ccode_fe114+ ccode_fe115+ ccode_fe116+ ccode_fe117+ ccode_fe118+ ccode_fe119+ ccode_fe120+ ccode_fe121+ 
ccode_fe122+ ccode_fe123+ ccode_fe124+ ccode_fe125+ ccode_fe126+ ccode_fe127+ ccode_fe128+ ccode_fe129+ 
ccode_fe130+ ccode_fe131+ ccode_fe132+ ccode_fe133+ ccode_fe134+ ccode_fe135+ ccode_fe136+ ccode_fe137+ 
ccode_fe138+ y_fe1+ y_fe2+ y_fe3+ y_fe5+ y_fe6+ 
y_fe7+ y_fe8+ y_fe9+ y_fe10+ y_fe11+ y_fe12+ y_fe13+ y_fe14+ y_fe15+ y_fe16+ y_fe17+ y_fe18+ y_fe19+ y_fe20+ 
y_fe21+ y_fe22+ y_fe23+ y_fe24+ y_fe25+ y_fe26+ y_fe27+ y_fe28+ y_fe29+ y_fe30+ y_fe31+ y_fe32+ y_fe33+ y_fe34+ 
y_fe35+ y_fe36+ y_fe37+ y_fe38+ y_fe39+ y_fe40+ y_fe41+ y_fe42+ y_fe43+ y_fe44+ y_fe45+ y_fe46+ y_fe47+ y_fe48+ 
y_fe49+ y_fe50+ y_fe51+ y_fe52+ y_fe53+ y_fe54+ y_fe55+ y_fe56+ y_fe57+ y_fe58+ y_fe59+ y_fe60+ y_fe61"

model_sar_boix <- milm(fml=FML_SAR_BOIX, midata=newdata)
model_sar_boix<- data.frame(model_sar_boix$terms, model_sar_boix$beta, model_sar_boix$SE, model_sar_boix$beta/model_sar_boix$SE)
colnames(model_sar_boix) <- c("Variable","Est. Coef.","SE","z-score")

print(paste("BMR SAR"))
print(model_sar_boix)


#Model BOIX - SLX
FML_SLX_BOIX <-"draw ~ draw.lag + SLreg_boix_t1 + l_boix +
l_gdp_gro + l_lngdppc + l_lnpop + civil_war + international_war + l_resdep + ccode_fe1+ 
ccode_fe2+ ccode_fe3+ ccode_fe4 + ccode_fe5+ ccode_fe6+ ccode_fe7+ ccode_fe8+ ccode_fe9+ 
ccode_fe10+ ccode_fe11+ ccode_fe12+ ccode_fe13+ ccode_fe14+ ccode_fe15+ ccode_fe16+ ccode_fe17+ 
ccode_fe18+ ccode_fe19+ ccode_fe20+ ccode_fe21+ ccode_fe22+ ccode_fe23+ ccode_fe24+ ccode_fe25+ 
ccode_fe26+ ccode_fe27+ ccode_fe28+ ccode_fe29+ ccode_fe30+ ccode_fe31+ ccode_fe32+ ccode_fe33+ 
ccode_fe34+ ccode_fe35+ ccode_fe36+ ccode_fe37+ ccode_fe38+ ccode_fe39+ ccode_fe40+ ccode_fe41+ 
ccode_fe42+ ccode_fe43+ ccode_fe44+ ccode_fe45+ ccode_fe46+ ccode_fe47+ ccode_fe48+ ccode_fe49+ 
ccode_fe50+ ccode_fe51+ ccode_fe52+ ccode_fe53+ ccode_fe54+ ccode_fe55+ ccode_fe56+ ccode_fe57+
ccode_fe58+ ccode_fe59+ ccode_fe60+ ccode_fe61+ ccode_fe62+ ccode_fe63+ ccode_fe64+ ccode_fe65+ 
ccode_fe66+ ccode_fe67+ ccode_fe68+ ccode_fe69+ ccode_fe70+ ccode_fe71+ ccode_fe72+ ccode_fe73+ 
ccode_fe74+ ccode_fe75+ ccode_fe76+ ccode_fe77+ ccode_fe78+ ccode_fe79+ ccode_fe80+ ccode_fe81+ 
ccode_fe82+ ccode_fe83+ ccode_fe85+ ccode_fe86+ ccode_fe87+ ccode_fe88+ ccode_fe89+ 
ccode_fe90+ ccode_fe91+ ccode_fe92+ ccode_fe93+ ccode_fe94+ ccode_fe95+ ccode_fe96+ ccode_fe97+ 
ccode_fe98+ ccode_fe99+ ccode_fe100+ ccode_fe101+ ccode_fe102+ ccode_fe103+ ccode_fe104+ ccode_fe105+ 
ccode_fe106+ ccode_fe107+ ccode_fe108+ ccode_fe109+ ccode_fe110+ ccode_fe111+ ccode_fe112+ ccode_fe113+ 
ccode_fe114+ ccode_fe115+ ccode_fe116+ ccode_fe117+ ccode_fe118+ ccode_fe119+ ccode_fe120+ ccode_fe121+ 
ccode_fe122+ ccode_fe123+ ccode_fe124+ ccode_fe125+ ccode_fe126+ ccode_fe127+ ccode_fe128+ ccode_fe129+ 
ccode_fe130+ ccode_fe131+ ccode_fe132+ ccode_fe133+ ccode_fe134+ ccode_fe135+ ccode_fe136+ ccode_fe137+ 
ccode_fe138+ y_fe1+ y_fe2+ y_fe3+ y_fe5+ y_fe6+ 
y_fe7+ y_fe8+ y_fe9+ y_fe10+ y_fe11+ y_fe12+ y_fe13+ y_fe14+ y_fe15+ y_fe16+ y_fe17+ y_fe18+ y_fe19+ y_fe20+ 
y_fe21+ y_fe22+ y_fe23+ y_fe24+ y_fe25+ y_fe26+ y_fe27+ y_fe28+ y_fe29+ y_fe30+ y_fe31+ y_fe32+ y_fe33+ y_fe34+ 
y_fe35+ y_fe36+ y_fe37+ y_fe38+ y_fe39+ y_fe40+ y_fe41+ y_fe42+ y_fe43+ y_fe44+ y_fe45+ y_fe46+ y_fe47+ y_fe48+ 
y_fe49+ y_fe50+ y_fe51+ y_fe52+ y_fe53+ y_fe54+ y_fe55+ y_fe56+ y_fe57+ y_fe58+ y_fe59+ y_fe60+ y_fe61"

model_slx_boix <- milm(fml=FML_SLX_BOIX, midata=newdata)
model_slx_boix<- data.frame(model_slx_boix$terms, model_slx_boix$beta, model_slx_boix$SE, model_slx_boix$beta/model_slx_boix$SE)
colnames(model_slx_boix) <- c("Variable","Est. Coef.","SE","z-score")

print(paste("BMR SLX"))
print(model_slx_boix)


#Model BOIX - Spatial Durbin Model (SDM)
FML_SDM_BOIX <-"draw ~ draw.lag + draw.lag.SLfariss + SLreg_boix_t1 + l_boix +
l_gdp_gro + l_lngdppc + l_lnpop + civil_war + international_war + l_resdep + ccode_fe1+ 
ccode_fe2+ ccode_fe3+ ccode_fe4 + ccode_fe5+ ccode_fe6+ ccode_fe7+ ccode_fe8+ ccode_fe9+ 
ccode_fe10+ ccode_fe11+ ccode_fe12+ ccode_fe13+ ccode_fe14+ ccode_fe15+ ccode_fe16+ ccode_fe17+ 
ccode_fe18+ ccode_fe19+ ccode_fe20+ ccode_fe21+ ccode_fe22+ ccode_fe23+ ccode_fe24+ ccode_fe25+ 
ccode_fe26+ ccode_fe27+ ccode_fe28+ ccode_fe29+ ccode_fe30+ ccode_fe31+ ccode_fe32+ ccode_fe33+ 
ccode_fe34+ ccode_fe35+ ccode_fe36+ ccode_fe37+ ccode_fe38+ ccode_fe39+ ccode_fe40+ ccode_fe41+ 
ccode_fe42+ ccode_fe43+ ccode_fe44+ ccode_fe45+ ccode_fe46+ ccode_fe47+ ccode_fe48+ ccode_fe49+ 
ccode_fe50+ ccode_fe51+ ccode_fe52+ ccode_fe53+ ccode_fe54+ ccode_fe55+ ccode_fe56+ ccode_fe57+
ccode_fe58+ ccode_fe59+ ccode_fe60+ ccode_fe61+ ccode_fe62+ ccode_fe63+ ccode_fe64+ ccode_fe65+ 
ccode_fe66+ ccode_fe67+ ccode_fe68+ ccode_fe69+ ccode_fe70+ ccode_fe71+ ccode_fe72+ ccode_fe73+ 
ccode_fe74+ ccode_fe75+ ccode_fe76+ ccode_fe77+ ccode_fe78+ ccode_fe79+ ccode_fe80+ ccode_fe81+ 
ccode_fe82+ ccode_fe83+ ccode_fe85+ ccode_fe86+ ccode_fe87+ ccode_fe88+ ccode_fe89+ 
ccode_fe90+ ccode_fe91+ ccode_fe92+ ccode_fe93+ ccode_fe94+ ccode_fe95+ ccode_fe96+ ccode_fe97+ 
ccode_fe98+ ccode_fe99+ ccode_fe100+ ccode_fe101+ ccode_fe102+ ccode_fe103+ ccode_fe104+ ccode_fe105+ 
ccode_fe106+ ccode_fe107+ ccode_fe108+ ccode_fe109+ ccode_fe110+ ccode_fe111+ ccode_fe112+ ccode_fe113+ 
ccode_fe114+ ccode_fe115+ ccode_fe116+ ccode_fe117+ ccode_fe118+ ccode_fe119+ ccode_fe120+ ccode_fe121+ 
ccode_fe122+ ccode_fe123+ ccode_fe124+ ccode_fe125+ ccode_fe126+ ccode_fe127+ ccode_fe128+ ccode_fe129+ 
ccode_fe130+ ccode_fe131+ ccode_fe132+ ccode_fe133+ ccode_fe134+ ccode_fe135+ ccode_fe136+ ccode_fe137+ 
ccode_fe138+ y_fe1+ y_fe2+ y_fe3+ y_fe5+ y_fe6+ 
y_fe7+ y_fe8+ y_fe9+ y_fe10+ y_fe11+ y_fe12+ y_fe13+ y_fe14+ y_fe15+ y_fe16+ y_fe17+ y_fe18+ y_fe19+ y_fe20+ 
y_fe21+ y_fe22+ y_fe23+ y_fe24+ y_fe25+ y_fe26+ y_fe27+ y_fe28+ y_fe29+ y_fe30+ y_fe31+ y_fe32+ y_fe33+ y_fe34+ 
y_fe35+ y_fe36+ y_fe37+ y_fe38+ y_fe39+ y_fe40+ y_fe41+ y_fe42+ y_fe43+ y_fe44+ y_fe45+ y_fe46+ y_fe47+ y_fe48+ 
y_fe49+ y_fe50+ y_fe51+ y_fe52+ y_fe53+ y_fe54+ y_fe55+ y_fe56+ y_fe57+ y_fe58+ y_fe59+ y_fe60+ y_fe61"

model_sdm_boix <- milm(fml=FML_SDM_BOIX, midata=newdata)
model_sdm_boix<- data.frame(model_sdm_boix$terms, model_sdm_boix$beta, model_sdm_boix$SE, model_sdm_boix$beta/model_sdm_boix$SE)
colnames(model_sdm_boix) <- c("Variable","Est. Coef.","SE","z-score")

print(paste("BMR SDM"))
print(model_sdm_boix)